Formation of Magnetic Microphases in Ca 3 Co 2 6 
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We study a frustrated quantum Ising model relevant for Ca 3 Co 2 6 that consists of a triangular lattice of 
weakly-coupled ferromagnetic (FM) chains. According to our quantum Monte Carlo (QMC) simulations, the 
chains become FM and form a three- sublattice "up-up-down" structure for T < T C i. In contrast, long-period 
spin-density- wave (SDW) microphases are stabilized along the chains for T C \ < T < T c . Our mean field 
solutions reveal a quasi-continuous temperature dependence of the SDW wavelength, implying the existence of 
metastable states that explain the very slow dynamics observed in Ca 3 Co 2 6 . We also discuss implications of 
microphases for the related multiferroic compounds Ca 3 CoMn0 6 and Lu 2 MnCo0 6 . 

PACS numbers: 75.10.Jm 75.25.-j 75.40.Mg 



Introduction. — Geometric frustration, low-dimensionality, 
and quantum fluctuations can lead to exotic phase transitions 
and states of matter d 13 such as the field-induced magne- 
tization plateaus of SrCu2(B03)2 GHS), the spin-driven "ne- 
matic" transition in pnictides (6H§), and dimensional reduc- 
tion in BaCuSi206 ||T0j |TT]|. Ca3Co206 is another example 
comprising a triangular lattice of FM Ising chains coupled by 
weak antiferromagnetic (AFM) exchanges. This compound 
exhibits field-induced magnetization steps whose heights de- 
pend on the field sweep history and rate lfT2lfT6l . We will 
show that this out-of-equilibrium behavior has its roots in ex- 
otic equilibrium properties that can be extended to the related 
multiferroic compound Ca3Co2- x Mn x 06 lfT7H22li . 

The Co 3+ ions (Co II) on the trigonal prism sites of 
Ca3Co206 contain 3d 6 localized electrons that generate an 
S = 2 spin with large Ising-like anisotropy | 23 -26|. These 
ions form a triangular lattice of FM Ising chains along the c- 
axis (Fig. [T} and the structure comprises three sublattices of 
layers stacked along the c-axis in an ABC ABC . . . configura- 
tion. Although the AFM inter-chain couplings J 2 and J3 |27) 
[Fig. [TJa)] are an order of magnitude smaller than the intra- 
chain FM exchange, = 2 x 10 K (23l|27l, we will show 
that they strongly affect the intra-chain spin correlations over 
a window of temperatures below T c . 

The initial interest in Ca3Co206 was triggered by the ob- 
servation of out-of-equilibrium magnetization steps measured 
below ~ 8K and ~ 3.6 T that appear at regular field inter- 
vals. Previous works invoked a "rigid-chain model": ev- 
ery chain is replaced by a single Ising spin by assuming 
T <sc \Ji\ [28-32]. Each spin of the resulting triangular lat- 
tice Ising model (TLIM), represents the magnetization of the 
whole chain and it is flipped if gid B H overcomes its molecular 
field. Within this simplified framework, the regular field in- 
tervals result from the equally- spaced discrete molecular field 
spectrum [28]. However, this 2D scenario was challenged 
by the recent discovery of long- wavelength intra-chain spin- 
density- wave (SDW) ordering below T c * 25 K (T6ll33ll34l 
Motivated by this discovery, Chapon initiated the study of a 
more realistic 3D lattice model by using a random-phase ap- 
proximation (RPA) which is only valid close to T = T c (351 . 

By combining QMC simulations and mean field (MF) solu- 



tions of the 3D quantum Ising model relevant for Ca3Co2C>6, 
we reproduce most of the measured zero-field properties. A 
sequence of soliton lattices that lead to the observed SDW 
order appears for Tqi < T < T c through the competi- 
tion between intra- and inter-chain couplings. While the 
transverse field stabilizes a ferrimagnetic (FIM) up-up-down 
(UUD) state below Tqi via order-by-disorder (36), very small 
longer-range exchange couplings, not included in our model, 
should be responsible for the actual T = ordering of 
Ca3Co2C>6 (37). Our MF solutions show that the ordering 
wave- vector changes quasi-continuously as a function of T, 
implying the existence of many competing metastable states. 
Even though the modulation wavelength increases for lower 
T and the rigid-chain picture is apparently applicable for 
T < Tqi, the relaxation is known to be extremely slow and 
practically never complete (T6l . According to our results, the 
observed slow dynamics for T < 10 K (T6l is a direct con- 
sequence of the multiple SDW microphases that appear for 
Tqi < T < T c . This exotic regime can only be captured by 
solving the 3D model beyond the RPA (35). 

Model — We use a pseudospin-1/2 to represent the lowest 
energy doublet (S z = ±2) of the Co II ions. The Hamiltonian 




FIG. 1. (Color online) (a) Exchange couplings between the Co II 
ions: The FM coupling J\ (a thick solid line) and the AFM couplings 
J 2 (thin solid lines) and J3 (dash lines). The lines within layers are 
projections of the inter-chain couplings, (b) The lattice projected on 
the ab plane. Each dot represents a chain. 




FIG. 2. (Color online) Upper panels show the specific heat while 
lower panels show |cr^ =0 ^ =0 | 2 for J 2 = J3 = 0.1|/i| and T = 0.3|/i|. 
OBC (PBC) is imposed in the c direction in (a) and (c) [(b) and (d)]. 
The arrows indicate a tiny anomaly in the low- T regime. The insets 
provide enlarged views. 
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where o~i is the vector of Pauli matrices for the ion i on the 
3D lattice shown in Fig. [T] J t j = J\ for nearest-neighbor 
(NN) sites along chains, and = J 2 (J3) for NN sites on 
NN (next-NN) layers [Fig.[TJa)]. Previous measurements in- 
dicate that J 1 is FM (Ji < 0) while J 2 , J3 \J\\ are AFM 
lfT4l . H = gciA B SB is the external magnetic field along the c- 
axis (g c is the gyromagnetic factor), while the transverse field 
is included for accelerating the QMC relaxation and remov- 
ing the macroscopic ground state degeneracy without invok- 
ing smaller and unknown longer-range exchange couplings. 
We will assume J 2 = J 3 = 0.1|/i|, T = 0.3|/i|, and H = 
unless otherwise specified (38). 

We use the continuous-time QMC method (39l|40| to com- 
pute the thermodynamic phase diagram. There is no sign 
problem because the frustrated terms are diagonal. The clus- 
ters can expand along each chain and the imaginary-time (r) 
direction. The weight factors due to J 2 and J3 appear in the 
cluster-flip attempt. We use the replica exchange method PHI 
for the lowest- T simulations l42l . The simulated lattice has 
Lx Lx L c unit cells with L c = 10L ( Mayer = 3L C layers). 

Figures [2ja) and [2jb) show the specific heat C obtained 
from our QMC simulations for periodic (PBC) and open 
boundary conditions (OBC) along the c-axis, while PBC are 
applied in the a and b directions. There are three different 
regimes. The A peak at T c 1.4|/i | indicates a 3D phase tran- 
sition from a paramagnetic phase to an ordered state. In ad- 
dition, there are two different ordered regimes below T c sep- 
arated by a tiny peak at T < 0.1|/i |, which we will refer to 
as intermediate- and low-7" regimes. While the position of 
this peak exhibits moderate size dependence, the consistent 
shift towards higher T for larger values of L implies robust- 



ness of the lowest- T phase against size effects. The sensitiv- 
ity of C(T) to the boundary conditions along the c-axis for 
T < T c is caused by a mismatch between the wave-vectors of 
the finite size lattice and the optimal SDW wave-vector in the 
intermediate- T regime. In what follows we adopt OBC along 
the c-axis because it is more convenient for detecting modula- 
tions with wavelength comparable to Mayer ( see below). 

FIM state in the low-T regime. — We will first discuss the 
state of equilibrium in the low-T regime. If T = T = 0, ev- 
ery chain is FM and the ground state subspace has the well- 
known massive degeneracy of the TLIM |43l . The lowest or- 
der correction to the ground state energy is 0(T 2 ). We intro- 
duce the hexagonal plaquette variables r^. = (1/2) Z</ 7 >^ cr z . 
with (ij)^ denoting sites connected by (ju = 2, 3). Since 

and 



any unperturbed ground state satisfies 0"?t^. 



(2)i ~ T \?>)i = T \ ' me ener gy cost °f flipping a spin of the ion i 
is AEi = 4|/i| + 4/|Tf I, where / = (J 2 + J 3 )/2. Therefore, the 
leading non-trivial contribution to the second-order effective 
Hamiltonian, <Hf = - 2,- T 2 /AE h is 
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Here we have used that the projection of |r^| in the 
unperturbed ground state subspace is a constant. Conse- 
quently, the lowest-order non-trivial effective interaction is 
a FM coupling between the next-NN chains that stabilizes 
the three- sublattice UUD state of FM chains. This FIM 
state has a spontaneous magnetization at 1/3 of the saturation 
value. To verify this numerically, we calculate (|cr^ 



|2\ - 



4=0,^=0' 



> = 



(lAT 1 ^ 1 2,- Jq 'dTof(T)| 2 >. Figures 2[c) and (d) show that 

d^l^o co=o\ 2 } a PP roacnes (1/3) 2 = 1/9 in the low-7 regime in 
agreement with our analytical result. The corresponding or- 
dering temperature coincides with the tiny anomaly in C(T). 

SDW order in the intermediate -T regime. — We will next 
discuss the most important intermediate- T regime. Fig- 
ure [3] shows the equal-time structure factor S(q) = 

r^n-^'"^^) slightly below T c (at T = 1.3|/i|), 

j 1 j 

which was extracted from the major peak of C(T). The Bragg 
peak at q 3 = Q3 is slightly shifted from q 3 = 2n/3 indi- 
cating modulated spin ordering along the c-axis [Fig. [3] (a)], 
while Fig. |3jb) clearly shows that each layer is FM. This 
is a three- sublattice SDW order with a relative phase shift 
of 27r/3; the numerical results are consistent with (o~ z iv ) ~ 
acos(q'r it 3 + <p v ) (v e {A,B,C}) in the single-harmonic ap- 
proximation, where q' = Q3 - 2n/3, r^ is the layer index 
of site i, and (pc - (pB = <Pb - <Pa = 2n/3. The very small 
value of \q'\ implies that the modulation period is very long: 
27r|g / |~ 1 - 3 x 10 2 at T = 1.3|/i|. These features become even 
more evident in the correlation functions of average moments 
per layer m/ = L" 2 2z;r /3 =/ (Fig.|4|). The abrupt decay close 
to the edges is a consequence of OBC. Although larger sys- 
tems are necessary to determine the precise T-dependence of 
q\ the obtained long-period modulation is in excellent agree- 
ment with recent experiments 1 16, 33ll34l. 
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FIG. 3. (Color online) S(q) at T = \3\Jx\ for J 2 = 7 3 = 0.1|/i|, 
T = 0.3|/i|, L = 16 (M ay e r = 480) and OBC along the c-axis. The 
vertical lines in (a) indicate q 3 = 2^/3 and 4tt/3, and the inset shows 
an enlarged view around q 3 = 2tt/3. The wave- vector is varied as (a) 
q = (0, 0, q 3 ) and (b) q = (qi,0,Q 3 ). Error bars are smaller than the 
symbol size. The line is a guide to the eye. 
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FIG. 4. (Color online) (m m/) for the same J 2 , J3, T, T, and L as in 
Fig. [3] OBC are imposed in the c direction and the center layer of the 
simulated lattice is chosen as / = 0. 



We now discuss the origin of the SDW and note that this or- 
dering does not appear in apparently similar lattices. For ex- 
ample, hexagonal lattice Ising systems, such as CSC0CI3 and 
CsCoBr 3 , exhibit the partially-disordered AFM state for inter- 
mediate T, even though they are also realizations of weakly- 
coupled Ising FM chains that form a triangular lattice l44ll . 
The crucial difference is in the connectivity of the inter-chain 
couplings: while the hexagonal lattice contains frustrated 
loops only within each layer, J 2 and J3 connect spins on dif- 
ferent layers (Fig. [TJ and consequently compete against the 
dominant intra-chain coupling J\ . 

Our numerical results suggest a natural MF approximation. 
We assume that each layer is FM as it is indicated by our QMC 
simulation. Indeed, the intra-layer effective FM coupling is 
induced not only by T but also by thermal fluctuations (it ap- 
pears in the second-order contribution of a high-7 expansion). 
In fact, our MC simulations show that the microphases still 
exist over an extended window of temperatures even in ab- 
sence of the transverse field. As expected, this phenomenon is 
entirely driven by the classical exchange interaction between 



Ising variables on a particular type of geometrically frustrated 
lattice. Therefore, we will take T = in what follows. The 
MF equations for the magnetization of each layer are 



(mi) = tanhyS/z/ 



(3) 



with hi = -/i«m/ +3 > + <m/_ 3 » - 3J 2 ((mi +1 ) + <m/_i» - 
3/3«m/+2) + (m/_2». The wave- vector Q c of the highest- T or- 
dered phase is given by the minimum of /mf(^) = 2/i cos 3q+ 
6J2 cos q + 6/3 cos 2q, and T c = -Jmf(Qc) E3. If J 2 and 
J3 are AFM, any finite inter-chain coupling leads to incom- 
mensurate SDW ordering at T = T c [35]. The MF solution 
for T < T c is obtained by solving Eq. ^ numerically. Fig- 
ure [5] shows q'(T) = Q(T) - 2n/3 (Q is the ordering vector) 
obtained by imposing PBC and by varying Mayer up to 2000 
{J 2 - h - 0.1|/i|). The q'(T) curve is qualitatively similar to 
the ones obtained for the ANNNI model |45|. However, we 
are not aware of any unambiguous realization of this prototyp- 
ical model in Mott insulators. The q f = phase (FM chains) 
stabilized at the lowest temperatures is the UUD FIM state ob- 
tained from our analytical approach and from the QMC simu- 
lations. This state becomes unstable for T > T C i ~ 2.171|/i | 
(the overestimation of Tqi is expected for a MF approxima- 
tion), above which q'(T) changes quasi-continuously. The 
obtained amplitude \q'\ is small in the entire regime, in agree- 
ment with experiments [16, 33 , 34 ] and with our QMC results. 
The optimal states are determined through close competition 
among many metastable states. This also implies that the fine 
structure of q'(T) should be very sensitive to small additional 
couplings that are not included in our model. However, the 
quasi-continuous change of q'(T) is a robust feature. 

A continuum approximation of our MF theory (analogous 
to Ref. 46) shows that the SDW phase corresponds to a quasi- 
continuous sequence of microphases driven by entropic ef- 
fects stabilizing a finite concentration of solitons (kinks) along 
the chains. The solitons form domain walls perpendicular to 
the c-axis. They crystallize into a lattice and the mean separa- 
tion A between walls determines q' oc The value of X is 
controlled by a balance between the chemical potential of soli- 
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FIG. 5. (Color online) q'{T) obtained from the MF theory (J 2 = J 3 
0.1 |). The inset shows an enlarged view. 
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FIG. 6. (Color online) Magnetization curve obtained by using a 
relaxation process described in the text, for J 2 = J3 = 0.1|/i|, 
r = 0.3|/i|, L = 12 (M ayer = 360), and OBC in the c direction. 

tons and an effective repulsive interaction that decays expo- 
nentially in the distance between solitons [46 1 . The outcome 
of this balance is that A diverges logarithmically in T - T C \ 
[our result (Fig. [5]) reproduces this behavior]. A number of 
metastable states appear in this regime with different mod- 
ulation periods. They are separated by free energy barriers 
associated with the creation/annihilation and redistribution of 
magnetic domain walls. These barriers give the dominant con- 
tribution to the observed slow dynamics because the relevant 
relaxation modes are suppressed at low T. 

Magnetization curve. — Finally, we present M(H) obtained 
in a simulated relaxation process in the realistic 3D model. 
We equilibrate the system at a given T for H = and then 
increase H gradually. We take 10 4 steps at each value of H, 
which is insufficient for equilibration at T «: T c . After reach- 
ing a sufficiently high field, we go back to H - in the same 
way and stop. For H ± 0, we only allow clusters to expand 
in the r direction, which corresponds to a classical single spin 
flip when T = 0. Although our dynamics is different from 
the real dynamics, our results reproduce the main experimen- 
tal observations, except for the less clear steps that appear 
at the highest-fields (above ~ 3.6 T) lfT2lfT6l . As is shown 
in Fig. ([6]), slightly below T c (T = 1.3|/i|), we only find a 
small feature suggesting a 1/3 plateau, which becomes more 
pronounced at T = 0.8|/i| accompanied by small hysteresis. 
Steps at regular magnetic field intervals appear at T = 0.3|/i|, 
which is still inside the SDW phase for H = 0. The reproduc- 
tion of equidistant steps in the relaxation dynamics supports 
the notion of metastability of the observed low-7 states. The 
heights of the steps obtained with the 3D model differ from 
the values obtained with the rigid-chain model I28H321 . 

Conclusions. — We reproduced the temperature dependent 
SDW state that was reported by recent neutron diffraction ex- 
periments in Ca3Co2C>6. More importantly, we showed that 
the SDW phase arises from a crystallization of domain walls 
that results in a large number of competing metastable states 
for Tqi < T < T c . By uncovering these microphases in 



Ca3Co206, we explained the origin of the extremely slow 
relaxation of the Bragg peaks fl6l . Disorder induced pin- 
ning of the domain walls that exist in the microphases also 
provides a natural explanation of the observed linear- T con- 
tribution in Ca3Co206 l47l . Order-by-disorder induced by 
a small transverse field leads to a FIM phase in the low-7 
regime. However, this result does not explain the recent obser- 
vation of an order-order transition to a different commensurate 
phase (371 . Therefore, although the FIM state is the ground 
state of Eq. ([T]), other subtle perturbations, such as intra-layer 
AFM exchange interactions between next-NN chains, must be 
included to explain the actual T - ordering of Ca3Co206. 

From our results we predict that microphases should also 
exist in the related multiferroic compounds Lu2MnCo06 l48l 
and Ca3Co2- x Mn x 06 (x « 1) |fT7H22l . Since magnetic do- 
main walls carry an internal electric dipole moment in these 
materials [49], the microphases should be sensitive to an ex- 
ternal electric field that introduces a bias between walls with 
opposite electric polarizations. Indeed, the dielectric constant 
of both compounds exhibits a broad peak below T c fT9ll48ll . 
We propose that this peak arises from the long-wavelength 
modulation of the electric dipole moments induced by differ- 
ent crystallization of magnetic domain walls (microphases). 

We thank S.-W. Cheong and T. Suzuki for valuable discus- 
sions. The numerical work was done in supercomputers of 
NERSC. Work at LANL was performed under the auspices 
of the U.S. DOE contract No. DE-AC52-06NA25396 through 
the LDRD program. 
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